Identifying hemodynamic factors associated with the rupture of anterior communicating artery aneurysms based on global modeling of blood flow in the cerebral artery network

Anterior communicating artery (ACoA) aneurysms are more prone to rupture compared to aneurysms present in other cerebral arteries. We hypothesize that systemic blood flow in the cerebral artery network plays an important role in shaping intra-aneurysmal hemodynamic environment thereby affecting the rupture risk of ACoA aneurysms. The majority of existing numerical studies in this field employed local modeling methods where the physical boundaries of a model are confined to the aneurysm region, which, though having the benefit of reducing computational cost, may compromise the physiological fidelity of numerical results due to insufficient account of systemic cerebral arterial hemodynamics. In the present study, we firstly carried out numerical experiments to address the difference between the outcomes of local and global modeling methods, demonstrating that local modeling confined to the aneurysm region results in inaccurate predictions of hemodynamic parameters compared with global modeling of the ACoA aneurysm as part of the cerebral artery network. Motivated by this finding, we built global hemodynamic models for 40 ACoA aneurysms (including 20 ruptured and 20 unruptured ones) based on medical image data. Statistical analysis of the computed hemodynamic data revealed that maximum wall shear stress (WSS), minimum WSS divergence, and maximum WSS gradient differed significantly between the ruptured and unruptured ACoA aneurysms. Optimal threshold values of high/low WSS metrics were determined through a series of statistical tests. In the meantime, some morphological parameters of aneurysms, such as large nonsphericity index, aspect ratio, and bottleneck factor, were found to be associated closely with aneurysm rupture. Furthermore, multivariate logistic regression analyses were performed to derive models combining hemodynamic and morphological parameters for discriminating the rupture status of aneurysms. The capability of the models in rupture status discrimination was high, with the area under the receiver operating characteristic curve reaching up to 0.9. The findings of the study suggest that global modeling of the cerebral artery network is essential for reliable quantification of hemodynamics in ACoA aneurysms, disturbed WSS and irregular aneurysm morphology are associated closely with aneurysm rupture, and multivariate models integrating hemodynamic and morphological parameters have high potential for assessing the rupture risk of ACoA aneurysms.


Introduction
Intracranial aneurysms (IAs), featured by abnormal local enlargement of cerebral arteries, are a common vascular disease affecting approximately 1%-5% of the population (Brisman et al., 2006;Starke et al., 2013).Rupture of an IA usually causes subarachnoid hemorrhage, leading to a mortality rate of up to 50% (Schröder et al., 1997;Van Gijn et al., 2007), with over 40% of survivors experiencing severe disabilities (Amenta et al., 2012).IAs can develop in multiple cerebral arteries, with about 15.5% located at the anterior communicating artery (ACoA) (Morita et al., 2012).Blood flow patterns in the ACoA are highly complex, because multiple arteries converge and diverge here, which renders the blood flow prone to disturbance and sensitive to flow conditions in adjacent arteries.The intensity of flow disturbance in the ACoA may be augmented when an aneurysm that significantly alters local vascular morphology develops, which may further increase the risk of endothelial dysfunction and vascular damage (Ku et al., 1985;Malek et al., 1999).In fact, it has been found that ACoA aneurysms are particularly more prone to rupture compared to aneurysms located elsewhere in patients with multiple IAs (Lu et al., 2013).Another critical issue associated with the management of ACoA aneurysms is the relatively high risk of complications or fatal events subsequent to endovascular treatments compared with IAs located elsewhere (Park et al., 2018;Chen et al., 2020).Therefore, identifying ACoA aneurysms of high rupture risk has great clinical significance.
Previous studies have demonstrated the potential associations of many morphological and hemodynamic factors (e.g., larger aneurysm size, large aspect ratio, abnormally low/high/oscillatory wall shear stress) with aneurysm rupture (Chien et al., 2008;Backes et al., 2014;Li et al., 2018;Liang et al., 2018).However, it remains debatable whether such associations could be applied to all IAs.In particular, it was found that many ACoA aneurysms ruptured even though their sizes were small (Forget et al., 2001;Jeong et al., 2009;Joo et al., 2009), which implies that large aneurysm size and size-dependent hemodynamic disturbances might not be reliable factors for predicting the rupture risk of ACoA aneurysms.In this context, some studies attempted to explore factors associated with the formation and rupture of ACoA aneurysm from a more systemic viewpoint, e.g., considering the anatomical structure of the circle of Willis (CoW)) and associated hemodynamic characteristics.A study found that unilateral hypoplasia of the A1 segment of the anterior cerebral artery (ACA1) was detected in 68.1% of patients with ACoA aneurysms (Kwak et al., 1980).In comparison, such anatomical defect was found in only 10%-35% of adults free of cerebral vascular disease (Niederberger et al., 2010).Some other studies revealed that anatomical asymmetries of the CoW were associated with ACoA aneurysm instability and rupture (Marc et al., 2012;Tarulli et al., 2014;Stojanović et al., 2019).In general, anatomical variations in CoW would not cause brain ischemia, but will alter the distribution of blood flow within the CoW and blood flow patterns in specific cerebral arteries (Alastruey et al., 2007;Liang et al., 2011).Blood flow in the ACoA or an aneurysm there would be particularly sensitive to CoW variations given the specific location of the ACoA in the CoW, as has been demonstrated by previous studies (Karmonik et al., 2009;Hassan et al., 2011;Liang et al., 2016).Based on these previous findings, we hypothesize that systemically considering the anatomical variations in CoW and their hemodynamic effects would favor a more sophisticated exploration of hemodynamic factors related to the rupture of ACoA aneurysms.
Medical image-based computational fluid dynamics (CFD) techniques have been widely applied to address hemodynamic differences between ruptured vs unruptured IAs thereby identifying hemodynamic factors of potential value for predicting the risk of aneurysm rupture (Cebral et al., 2011;Xiang et al., 2011;Li et al., 2014).Most studies in this field built CFD models covering a limited physical region (usually confined to an IA and its adjacent arteries).The modeling strategy can effectively lower the burden of model reconstruction and reduce computational cost, thereby facilitating studies on large-scale IAs.The advantages however may be compromised by the reduced fidelity of numerical results due to the missing information on systemic hemodynamics that is critical for determining the inflow and outflow conditions of IAs.At this point, some studies have proved the importance of including sufficient upstream cerebral arteries for accurate computation of intra-aneurysmal hemodynamics (Castro et al., 2006;Pereira et al., 2013;Hodis et al., 2015;Hua et al., 2015;Xin et al., 2022).A recent study on a small number of ACoA aneurysms demonstrated that the computed time-averaged wall shear stresses in ACoA aneurysms by unilateral model and complete model of the CoW could differ up to 62% in certain cases (Yang et al., 2024).
In light of the aforementioned issues, we firstly conducted numerical experiments to investigate to what extend the selection of physiological boundaries in the modeling of ACoA aneurysm would affect the computed intra-aneurysmal hemodynamic quantities (herein simplified models confined to the local aneurysm region vs global models of the entire CoW).As we will show later, the two modeling strategies resulted in differential numerical results.Subsequently, we applied the global modeling strategy to forty ACoA aneurysms (containing twenty ruptured and twenty unruptured ones), yielding a series of numerical data for evaluating the hemodynamic differences between ruptured and unruptured aneurysms.In addition to the improvement in modeling method, we detailed the analysis of hemodynamic data by expanding wall shear stress (WSS) metrics for comparison and refining the determinant of threshold values.Finally, we combined hemodynamic and morphological parameters to establish mathematical models for predicting the rupture risk of ACoA aneurysms.

Patients and clinical data
Computed tomography angiography (CTA) data were retrospectively collected from patients diagnosed with ACoA aneurysms at Renji Hospital (Shanghai, China).The use of the data for research has been approved by the hospital's ethics committee.In consideration of global modeling strategy adopted in our study, only the aneurysms with high-quality medical images supporting the reconstruction of all major cerebral arteries were included.As a consequence, forty ACoA aneurysms were involved in our study, containing twenty ruptured ones and twenty unruptured ones.Preliminary image analysis showed that hypoplastic (underdeveloped) ACA1 was present in seven patients in the ruptured aneurysm group and in six patients in the unruptured aneurysm group, which implies that hypoplasia of the ACA1 may not be a factor related to the rupture risk of ACoA aneurysm, at least in the cohort of patients involved in the present study.Boundary conditions assigned to the model of a patient (Pt1) (A), and the comparison of the simulated inflow waveforms for three patients (B).The inflow waveforms are those in the left ICA, right ICA, and BA.The RCR models are connected to the outlets of the six efferent arteries (i.e., left/right MCA, left/right ACA2, and left/right PCA2).Note that the added extension tubes are not shown in panel (A).

Image-based geometric model reconstruction and mesh generation
The geometric model of each ACoA aneurysm was reconstructed along with large cerebral arteries from CTA images using Mimics 16.0 (Materialise, Belgium).To control the complexity of model while maintaining the natural anatomical structure of the CoW, the following model reconstruction strategies were adopted: 1) all the large cerebral arteries constituting the CoW are maintained [i.e., A1 segment of left/right anterior cerebral artery (ACA1), anterior communicating artery (ACoA), P1 segment of left/right posterior cerebral artery (PCA1), and left/right posterior communicating artery (PCoA)]; 2) major afferent cerebral arteries [i.e., left/right internal carotid artery (ICA), basilar artery (BA)] are fully reconstructed based on available image data; and 3) reconstruction of efferent cerebral arteries were limited to the six largest ones [i.e., A2 segment of left/right anterior cerebral artery (ACA2), left/right middle cerebral artery (MCA), and P2 segment of left/right posterior cerebral artery (PCA2)] that deliver blood to different brain tissue regions.The reconstructed model for each patient was further treated by adding extension tubes to the inlets of the afferent arteries and the outlets of the efferent arteries.The length of the extension tube connected to each inlet/outlet was set to be 15 times/6 times of the vascular diameter (Myers et al., 2001).
The addition of the extension tubes was expected to minimize the influence of imposed pressure or flow velocity distribution at the inlets/outlets on the flow patterns in the aneurysm and CoW.
Mesh generation based on the reconstructed geometric model was implemented in Fluent Meshing (Fluent, ANSYS Inc., USA).Triangle surface meshes were firstly created on the vascular walls, followed by the division of the inner fluid domain with tetrahedral elements.Finally, ten layers of prism elements were mapped along the vascular wall to replace the tetrahedral elements in the near-wall regions.The thickness of the first prism layer adjacent to the wall was set to be 0.005 mm, which was increased by a ratio of 1.25 towards the inner layers.The prism element layers, given the regular mesh topology and smooth changes in mesh size from the wall towards the interior, were expected to improve the precision of flow simulation in the near-wall region where velocity gradient is large (Xu et al., 2018a).Mesh independence analysis was performed on a randomly selected model, by varying the maximum size of the tetrahedral elements from 0.16 mm to 0.14 mm and 0.12 mm in sequence.It is noted that the minimum-to-maximum mesh size ratio was fixed at 1/16.The results showed that the relative error of the computed time-averaged WSS (TAWSS) at the entire aneurysm wall decreased from 1.747% to 0.555% following the reduction in mesh size.Therefore, the maximum size of the tetrahedral elements was set at 0.14 mm and was applied to all the models.The resulting mesh models contained 8 million to 25 million elements, depending on the size and complexity of the geometric models.

Setup of computational model
Flow conditions imposed to the boundaries of a model have important influence on the outcome of numerical simulation (Xu et al., 2018b;Steinman and Pereira, 2019;Mehdi et al., 2021).However, blood flow velocities/rates in cerebral arteries are not routinely measured in clinical diagnosis or treatment of IAs (Xin et al., 2022).In this study, we employed a zero-one dimensional (0-1D) model of the cardiovascular system (Liang et al., 2011) to simulate the flow waveforms in cerebral arteries.The model was partly personalized by incorporating the detailed structural and geometrical information of the cerebral artery network (extracted from the 3D geometric model reconstructed based on CTA images), as such it could account for the influence of patient-specific cerebral artery anatomy while simulating flow waveforms in cerebral arteries.The simulated flow waveforms in the afferent cerebral arteries were then used to prescribe the inlet flow conditions of each 3D model.Figure 1A shows the model reconstructed for a patient with prescribed boundary conditions.The inflow waveforms simulated (with 0-1D models) for three patients are compared in Figure 1B.Considerable inter-patient differences in the magnitudes of the flow waveforms in the three afferent arteries were observed, reflecting the influence of patient-specific geometrical size and structure of the cerebral artery network.The simulated flow waveforms were similar in shape with in vivo measurements reported in previous studies (Reymond et al., 2009;Zhou et al., 2020).The outlets of the efferent cerebral arteries were each connected to a three-element Windkessel model (i.e., RCR model) which plays a role of controlling the flow division among the efferent cerebral arteries.The RCR model was formulated in form of Eq. 1, where Q and P represent the volume flow rate and pressure at the outlet, which are related via two resistances (R p , R d ) and a compliance (C) (Alastruey et al., 2008;Pahlevan et al., 2011).The total resistance [R T = R p + R d , with R p /R T = 0.2 (Stergiopulos et al., 1992)] in each outlet model was tuned so that the computed flow division among the efferent cerebral arteries is consistent with the populationaveraged physiological data (Liang et al., 2011).All vascular walls in the 3D model were assumed to be rigid to which the noslip boundary condition was imposed.
Blood flows in cerebral arteries and aneurysms were modeled as non-Newtonian incompressible fluid and laminar flow regime was assumed.Accordingly, blood flows were governed by the continuity equation (Eq.2) and momentum conservation equation (Eq.3).
where t is the time, u is the flow velocity vector, p is the pressure, and ρ (1,060 kg/m 3 ) is the blood density.The non-Newtonian rheological properties of blood were represented by modeling the dynamic viscosity as a function of shear rate (Eq.4, the Carreau model), where μ 0 (0.056 Pa•s) and μ ∞ (0.00345 Pa•s) stand for the dynamic viscosities of blood when the shear rate approaches zero and becomes infinitely large, respectively, λ and n are constants, with their values being set at 3.313 s and 0.3568, respectively (Johnston et al., 2004;Nagargoje et al., 2021).
The governing equations of blood flow were numerically solved along with the boundary conditions with the finite volume method using a CFD package (Fluent, ANSYS Inc., USA).Second-order schemes were adopted for both the temporal and spatial terms.The outlet models were solved in conjunction with the 3D model via userdefined functions (UDFs) embedded in Fluent.Specifically, the pressures at the 3D model outlets were updated via solving Eq. 1 (with the explicit Euler method) during internal iterations within each numerical time step.The numerical time step was fixed at 0.001s and the threshold of residuals for judging convergence at every time step was set to be 1E-4.Numerical simulations were conducted under pulsatile flow conditions.Each set of numerical simulation was run continually for two cardiac cycles to achieve periodic solutions.

Data analysis
Analysis of hemodynamic quantities was focused on WSS at aneurysm wall.Under pulsatile flow conditions, WSS vector (τ) changes in magnitude and direction.To facilitate quantitative analysis, various WSS metrics have been defined in the literature.Typical WSS metrics include the time-averaged WSS (TAWSS, averaged magnitude of τ over a cardiac cycle) and oscillatory shear index (OSI, a measure of directional changes in WSS over a cardiac cycle) (Xin et al., 2022).In the study, we further considered some metrics characterizing the spatial variations of WSS, since WSS topology has been demonstrated to affect the initiation and progression of vascular disease.For instance, negative WSS divergence, which is often detected in wall regions where the wall shear stress vectors converge toward a focal point, can predict the localization of atheroprone low-density lipoprotein (LDL) deposition (Mazzi et al., 2022;Tian et al., 2023), whereas high WSS gradient and gradient oscillation are associated with endothelium injury and aneurysm formation (Shimogonya et al., 2009;Metaxa et al., 2010).Herein, three spatial WSS metrics were considered, namely, WSS divergence (WSSD, Eq. 5), WSS gradient (WSSG, Eq. 6), and gradient oscillatory number (GON, Eq. 7), which were firstly calculated at each time step and then integrated over a cardiac cycle.
To further facilitate the comparison of WSS metrics between ruptured and unruptured aneurysms, the maximum, minimum and mean values of each WSS metric over the aneurysm wall were derived, which are denoted, for example, as maxTAWSS, minTAWSS, and SA-TAWSS for TAWSS.Furthermore, we calculated the total area of aneurysm walls exposed to WSS metrics higher or lower than certain thresholds.For TAWSS, AHS and ALS represent the areas of local aneurysm walls exposed to TAWSS > WSS HT or TAWSS < WSS LT , where WSS HT and WSS LT are the thresholds of high TAWSS and low TAWSS, respectively.Similar area indices were defined for high OSI (AHI, OSI > OSI HT ), low WSSD (ALD, WSSD < WSSD LT ), high WSSG (AHG, WSSG > WSSG HT ) and high GON (AHN, GON > GON HT ), respectively.Each area index was further divided by the total area of the aneurysm wall to derive the relative area index (i.e., area ratio).It is noted that the thresholds of each WSS metric were not derived directly from the literature or estimated empirically, but determined via a series of statistical tests so that the area indices differ most evidently between ruptured and unruptured aneurysms.
The Mann-Whitney U test was performed to compare parameters between ruptured and unruptured aneurysms.All the tests were two-tailed and a p value of <0.05 was used to judge

Influence on numerical results of physical boundary truncation strategies in image-based model reconstruction
We randomly selected two patients (Pt4 with a ruptured ACoA aneurysm and Pt5 with an unruptured ACoA aneurysm) for test.For each patient, two geometric models were reconstructed based on medical images with different physical boundary truncation strategies: one is a global model containing the ACoA aneurysm and all the major cerebral arteries, and the other one is a local model containing the ACoA aneurysm and the surrounding arteries (i.e., left and/or right ACA1, and left/right ACA2).The outlets of the two models were supported by the same RCR models, while the flow waveforms imposed at the inlets of the local model were derived from the numerical results of the global model.The boundary condition settings guaranteed that the global and local models had the same flow rates through the ACA1s and flow division between the left and right ACA2s.
Figure 2 shows the contour maps of WSS metrics in the aneurysm region and flow velocities in six representative vascular cross sections computed by the global and local models of a patient (Pt4).Visually, all the WSS metrics differed between the two models in terms of magnitude and spatial distribution.As for the flow velocity contours, evident between-model differences were observed for those in the cross sections of the left ACA1 and right ACA1, especially at the ACA1inlets (i.e., cross sections 'a-a' and 'b-b'), whereas the differences were relatively small in the cross sections of the left and right ACA2s.Mechanisms underlying the observations are related to the role of upstream arteries in shaping the 3D patterns of flows entering the ACA1s.In comparison with the global model, the local model, though imposed with the same volume flow rates in the ACA1s, lost the 3D flow information in the upstream arteries that would pose influence on flow patterns in the ACA1s, which ultimately resulted in differential predictions of WSS metrics in the ACoA aneurysm.To enhance quantitative comparison between the outcomes of the two types of models, the space-averaged values of some WSS metrics over the ACoA aneurysm wall were calculated and provided in Table 1.The between-model differences were within ±10% for SA-TAWSS and SA-WSSG, but were much larger for SA-WSSD.In addition, the between-model differences differed significantly between the two patients.These results indicate that errors in numerical results associated with local modeling are somehow uncertain and are dependent on patient-specific anatomical features of the cerebral artery network and ACoA aneurysm.In this sense, global modeling may be a reliable approach to minimizing the influence of improper physical boundary truncation, thereby improving the physiological fidelity of numerical results, although the computational cost will increase significantly compared with local modeling.The modeling strategy was adopted for all the ACoA aneurysms involved in the study.

Comparisons of morphological and hemodynamic parameters between ruptured and unruptured aneurysms
Hemodynamic simulations were performed for all the forty ACoA aneurysms.Figure 3 shows the contour maps of modelpredicted WSS metrics in four ACoA aneurysms (two unruptured vs two ruptured).The spatial distributions of the WSS metrics were highly complex, exhibiting aneurysm-specific characteristics.However, no specific characteristics were found to clearly differ between the unruptured and ruptured aneurysms.
Table 2 shows the statistical comparisons of the morphological and hemodynamic parameters between the twenty ruptured aneurysms and twenty unruptured aneurysms.The data are presented in form of mean ± standard deviation.Four out of nine morphological parameters were found to differ statistically between the ruptured and unruptured aneurysms, which included the nonsphericity index (NSI), height-to-width ratio (HWR), aspect ratio (AR), and bottleneck factor (BF).All the parameters were larger in the ruptured aneurysm group.Interestingly, the parameters measuring the absolute size of aneurysm, such as the dome area, volume, and aneurysm height, did not show statistically significant between-group differences.As for hemodynamic parameters, three WSS metrics (i.e., maxTAWSS, maxWSSG and minWSSD) were found to differ significantly between the unruptured and ruptured aneurysms.Specifically, maxTAWSS and maxWSSG were larger while minWSSD was smaller in the ruptured aneurysm group.The ruptured aneurysm group had an overall lower minTAWSS, but the between-group differences did not reach statistical significance.No between-group difference was identified for any space-averaged (SA) WSS metrics.
Although the space-averaged WSS metrics did not exhibit significant between-group differences, there still exists the possibility that WSS metrics limited to certain ranges differ between groups and could be used to discern the rupture status of ACOA aneurysms.Therefore, we further performed statistical analyses to explore the high or low thresholds of WSS metrics (herein, TAWSS, OSI, WSSD, WSSG, and GON).To facilitate quantitative analysis, the area and area ratio (relative to the wall area of aneurysm) of walls exposed to WSS metric values lower or higher than a given threshold (low threshold (LT), or high threshold (HT)) were calculated, which were subsequently compared between the ruptured and unruptured aneurysms by means of Mann-Whitney U test.To ensure that the thresholds are determined in a robust way, the candidate threshold of each WSS metric was varied in a large range, with the p value outputted from the Mann-Whitney U test being used as the criterion for judging the performance of the threshold.The p values are plotted against the varying thresholds in Figure 4.A valley was discernible on the p value-threshold curves of most WSS metrics, indicating the existence of an optimal threshold at which the p value reaches a local minimum.Notably, WSSD and WSSG exhibited a broad range of thresholds at which the corresponding p values are below 0.05, which implies that the power of ALD (area of low WSSD) and AHG (area of high WSSG) in discriminating the rupture status of ACoA aneurysms is less sensitive to the choice of threshold.In contrast, no thresholds resulting in a p value lower than 0.05 were found for low TAWSS, high OSI, and high GON.As for high TAWSS, the p value decreased to <0.05 when the threshold value was in the range of 8-10 Pa, but again increased over 0.05 following further increase of the threshold value.Based on the results, five area indices, namely, AHS (area of high TAWSS), ALD, AHG, ALDR (ALD ratio), and AHGR (AHG ratio), were found to differ significantly between the ruptured and unruptured aneurysms.The thresholds of the corresponding WSS metrics (i.e., TAWSS, WSSD, WSSG) were set at TAWSS HT = 8.3 Pa, WSSD LT = −6 Pa/m, and WSSG HT = 55,000 Pa/m, yielding p values of 0.043, 0.002, 0.003, 0.008, and 0.009, respectively.

Performances of individual morphological and hemodynamic parameters for discriminating the rupture status of aneurysm
Receiver operating characteristic (ROC) analysis was carried out to evaluate the capability of each individual morphological/hemodynamic parameter in discriminating the rupture status of aneurysms.Table 3 summarizes the results of ROC analyses for all the parameters.maxWSSG exhibited the highest capacity of discernment, achieving an area under the ROC curve (AUC) of 0.7925, followed by ALD, minWSSD, AHG, NSI, and AR whose AUCs were all above 0.75.If the sensitivity-specificity balance was further considered, maxWSSG, minWSSD, AHG, and NSI mildly outperformed ALD and AR.

Derivation of multivariate models for predicting the rupture status of aneurysm
We further explored whether combining different parameters could enhance the discrimination of aneurysm rupture status beyond individual parameters by means of multivariate logistic regression analysis.In order to minimize the influence of multicollinearity caused by variable intercorrelations on the outcome of regression analysis, Spearman's rank correlation analysis was firstly carried out to identify independent variables.The results are summarized in Table 4, where it was observed that maxWSSG and AHG were closely intercorrelated, and were correlated strongly with many other hemodynamic parameters, such as maxTAWSS, minWSSD, ALD, and ALDR.The two key morphological parameters (i.e., NSI and AR, which were identified by single-variate ROC analysis to be associated with aneurysm rupture) were independent, and were uncorrelated or correlated weakly with maxWSSG and AHG.Based on the above results, we selected two sets of parameters (one set consisting of maxWSSG, NSI and AR, the other set consisting of AHG, NSI and AR) as the independent variables in the multivariate logistic regression analysis.Given the twenty aneurysm rupture events in the dataset, the number of events per variable (EPV) (20 divided by 3) was about 7, which surpasses the lower limit of EPV (≥5) in logistic regression analysis for binary problem (herein, discriminating the ruptured status of aneurysm) (Peduzzi et al., 1996;Vittinghoff and McCulloch, 2007).To reduce potential roundoff errors caused by the differential magnitudes of parameters in multivariate regression analysis, each parameter was normalized by its mean value in all the aneurysms (herein the normalized parameter is denoted by a '*' subscript).In addition, since parameter values might be nonlinearly related to the rupture status of aneurysm, single-variate logistic regression analysis was conducted for each parameter to determine the order of polynomial regression model with the best predictive performance (evaluated with ROC analysis).The results indicated that the optimal orders of singlevariate polynomial models were second, fifth, fourth, and fifth for maxWSSG, AHG, NSI, and AR, respectively.However, when highorder polynomials of all the parameters are included in a multivariate model, the degrees of freedom (i.e., the number of predictors) will increase significantly, leading to low sample size corresponding to each individual predictor thereby augmenting the risk of overfitting (Babyak, 2004).With this in mind, we set the highest order of polynomials at second for maxWSSG and AHG, whereas first for NSI and AR in multivariate logistic regression analyses based on the trends of changes in AUC with the polynomial order in the single-variate logistic regression analyses so that the sample size was 10 [40 (number of aneurysms) divided by 4 (number of predictors)], fulfilling the lower limit of sample size (>8) in regression analysis as recommended by a recent study (Jenkins and Quintana-Ascencio, 2020).The multivariate logistic regression analyses with the two sets of parameters yielded two models predicting the probability of aneurysm rupture.
For maxWSSG, NSI and AR, the regression model was expressed as odd rp,1 exp −12.268 + 2.638maxWSSG* − 0.374 maxWSSG* Whereas for AHG, NSI and AR, the regression model was written as odd rp,2 exp −14.027 + 2.927AHG* − 0.854 AHG* ( ) 2 + 10.076NSI* + 3.524AR*] (9) By definition, the probability of aneurysm rupture was calculated as When the two multivariate regression models were examined with ROC analysis, the values of AUC were 0.9025 and 0.8725, respectively.The detailed results of ROC analysis for the multivariate regression models are listed together with those for the single-variate regression models in Table 5, and the corresponding ROC curves are plotted in Figure 5.It is clear that the multivariate models significantly outperformed any singlevariate models, which demonstrates the potential of combining multiple hemodynamic and morphological parameters through multivariate modeling for improving the prediction of aneurysm rupture.

Discussion
While numerous studies have been conducted to explore the associations between hemodynamic factors and the rupture status/ event of IAs, the outcomes of many studies might have been compromised by improper hemodynamic modeling.The issue could be especially important for IAs with multiple surrounding arteries, such as ACoA aneurysms.Theoretically, the blood flow patterns in an ACoA aneurysm are determined not solely by the aneurysm morphology but affected considerably by the multistream inflows and outflows as well.Our numerical tests on two ACoA aneurysms demonstrated that in comparison with global modeling where an ACoA aneurysm is modeled as part of the systemic cerebral artery network, local modeling of an ACoA aneurysm with the physical boundaries confined to its adjacent arteries could not accurately capture intra-aneurysmal hemodynamic characteristics even though the inflow rates were correctly assigned.Inspired by this finding, we carried out hemodynamic simulations for forty ACoA aneurysms (including ROC curves of different logistic regression models: maxWSSG, AHG, NSI, AR, maxWSSG-NSI-AR, AHG-NSI-AR. Frontiers in Bioengineering and Biotechnology frontiersin.orgtwenty ruptured and twenty unruptured aneurysms) based on the global modeling method.Although the global modeling method incurred remarkably increased computational cost, the numerical results would better reflect the in vivo hemodynamic conditions.Based on the numerical results, we explored hemodynamic parameters (herein mainly WSS metrics) associated closely with the rupture status of ACoA aneurysms, and subsequently combined the identified hemodynamic parameters with morphological parameters to derived models for predicting the rupture status of ACoA aneurysms via multivariate logistic regression analysis.
The maximum/minimum values (i.e., maxTAWSS, minWSSD, and maxWSSG) and area indices (AHS, ALD, AHG, ALDR, and AHGR) of TAWSS, WSSD, and WSSG were found to differ significantly between the ruptured and unruptured ACoA aneurysms.The finding regarding the close association of maxTAWSS with aneurysm rupture is consistent with previous findings (Detmer et al., 2019;Xin et al., 2022).Relatively, WSSD and WSSG were less addressed by previous studies.Our study demonstrated that negative WSSD and high WSSG were associated closely with aneurysm rupture as well.Physically, negative WSSD is detected when WSS vectors converge toward a focal region, which has been proved to enhance the accumulation of atherogenic substances (e.g., LDL) on vascular lumen surface (Mazzi et al., 2022;Tian et al., 2023).High WSSG indicates the presence of large spatial variations in the magnitudes of WSS vectors, which has been found to have the effect of upregulating inflammatory marker expression (Rouleau et al., 2010).Given the negative impact of negative WSSD and high WSSG on mass transfer and vascular function, aneurysm walls exposed to negative WSSD and/or high WSSG may be prone to functional and structural degenerations that might increase the risk of rupture.On the other hand, correlation analyses revealed that most of the WSS metrics identified to be associated with aneurysm rupture were strongly intercorrelated, although they have differential physical meanings, which implies that ruptured ACoA aneurysms might suffer from multiple coexisting patterns of WSS disturbance.Interestingly, low WSS and high OSI were not found to be associated with aneurysm rupture, although they have been widely employed to predict the rupture risk of IAs (Liang et al., 2018), including ACoA aneurysms (Detmer et al., 2019).Reasons underlying the discrepancy between our findings and previous ones remain unknown and may involve multiple factors, such as the location of aneurysm, composition of samples, and methods of hemodynamic modeling, which are usually study-specific.In addition, the relatively small size of ACoA aneurysms at rupture compared to other aneurysms (Forget et al., 2001;Jeong et al., 2009;Joo et al., 2009) may be a contributing factor, since small bifurcation aneurysms provide limited sac space for the development of flow recirculation and separation to generate low and oscillatory WSS.
As for morphological factors, four morphological parameters (i.e., NSI, HWR, AR, and BF) were identified to be associated closely with aneurysm rupture, among which NSI and AR exhibited the strongest association.The results were basically consistent with the findings of previous studies (Liang et al., 2018;Detmer et al., 2019;Xin et al., 2022), suggesting that irregular shape predisposes aneurysms to rupture.An exception was the volume of aneurysm (V aneu ), which did not exhibit a statistically significant difference between the ruptured and unruptured ACoA aneurysms, although it has been suggested to be a general risk factor for aneurysm rupture (Detmer et al., 2019;Xin et al., 2022).A potential explanation for this disparity may be similar to that for the irrelevance of ACoA aneurysm rupture to low WSS and high OSI.The specific location of ACoA aneurysm and the small size at the moment of rupture determine that the inflow conditions and irregular aneurysm shape (e.g., measured by NSI and AR), rather than the size of aneurysm, dominate the occurrence of aneurysm rupturepromoting flow disturbances.
Most of the identified hemodynamic and morphological parameters could independently discriminate the rupture status of aneurysms with an AUC of > 0.7 in ROC analysis, and the predictive capability of each individual parameter could be further improved (with the AUC increasing up to 0.85) if a nonlinear model was constructed via single-variate logistic regression analysis.Nevertheless, synthesizing multiple independent hemodynamic and morphological parameters in a single model by multivariate logistic regression analysis could yield better outcomes of prediction, e.g., AUC = 0.9025 for the combination of maxWSSG, NSI and AR, and AUC = 0.8725 for the combination of AHG, NSI and AR.These results suggest that hemodynamic and morphological factors are complementary in assessing the rupture risk of ACoA aneurysms.In addition, although the maxWSSG-NSI-AR model had larger AUC, we prefer to recommend the use of the AHG-NSI-AR model because maxWSSG is a local value susceptible to the influence of errors in model reconstruction and boundary condition prescription, whereas, AHG, which reflects the total area of aneurysm walls exposed to high WSSG, is relatively insensitive to such errors.
While the results of the present study provide valuable insights for exploiting the hemodynamic and morphological characteristics to assess the rupture risk of ACoA aneurysms, they must be considered in the context of certain limitations.A major limitation is due to the small data size.Although the forty ACoA aneurysms involved in the study enabled statistical analysis, the number of aneurysms could not support the derivation and testing of higher-order multivariate regression models for predicting rupture events given the requirement on sample size in regression analysis.Ideally, two datasets will be required if largerscale data would be available, with one dataset for deriving while the other dataset for testing models.From the results of single-variate regression analysis, higher-order models exhibited better performances than the first-order models, in this sense, there exists the potential for further improving the performances of multivariate regression models by including more high-order terms, provided sufficient data would be available.Another limitation is related to the rigid-wall assumption.The walls of cerebral arteries and IAs are compliant in nature, which deform dynamically under pulsatile flow conditions.Previous studies comparing rigid and compliant IA models have demonstrated that rigid models may considerably underestimate near-wall flow velocity oscillation and overestimate TAWSS in some regions (Bazilevs et al., 2010;Torii et al., 2011;Tupin et al., 2020;Sun et al., 2022).Fortunately, the discrepancies between the two types of models were small (usually <10%) in terms of the predicted results of key WSS metrics (e.g., space-averaged TAWSS, maxTAWSS).In consideration of the overall small size of ACoA aneurysms investigated by the present study, the impact of wall deformation on intra-aneurysmal hemodynamic characteristics should be secondary and would not pose substantial influence on the major findings of our study.Nonetheless, accounting for the influence of wall deformation with advanced modeling methods, such as the fluid-structure interaction method, would be expected to further improve the physiological fidelity of numerical results.In addition, we could not identify the specific sites of ruptured aneurysm walls based on CTA images and therefore could not explore the specific hemodynamic features at the rupture sites.This limitation may be addressed if other medical images (e.g., images taken by intracranial photography during surgery (Li et al., 2018) that enable the identification of rupture site were available).

Conclusion
A computational model-based study has been carried out to investigate hemodynamic differences between ruptured and unruptured ACoA aneurysms, and derive multivariate models for discriminating the rupture status of aneurysms.The study firstly demonstrated the importance of modeling ACoA aneurysm as part of the cerebral artery network in order to accurately capture intraaneurysmal hemodynamic characteristics, offering a methodological reference for relevant studies in this field.Indices of WSSD and WSSG, rather than traditional WSS metrics like low WSS area and high OSI, were found to differ significantly between the ruptured and unruptured ACoA aneurysms, which implies that ACoA aneurysms, which often rupture at small size, may differ from aneurysms in other cerebral arteries with respect to flow disturbances driving aneurysm wall weakening and rupture.In addition, multivariate regression models composed of polynomials of hemodynamic and morphological parameters were proved to outperform single-variate models in discriminating the rupture status of ACoA aneurysms, suggesting that the rupture risk of ACoA aneurysms should be assessed on the basis of a synthetical consideration of hemodynamic and morphological factors.On the other hand, it should be stressed that the derivation of multivariate regression models was limited by the small data size (herein n = 40), which could be further improved if larger-scale data were available.

FIGURE 2
FIGURE 2 Contour maps of WSS metrics in the aneurysm region and flow velocities in six representative vascular cross sections computed by a global model (A) compared with those computed by a local model (B).The models are constructed for a randomly selected patient (Pt4).The flow velocity contours are drawn based on the time-averaged flow velocities over a cardiac cycle.It is noted that the cross-sectional shapes of the ACA1 inlets differ slightly between the local and global models because the ACA1 inlet surfaces were slightly smoothened to facilitate the connection to the extension tubes in the case of local modeling.

FIGURE 4 p
FIGURE 4 p values of Mann-Whitney U test for comparing the ruptured and unruptured aneurysms in terms of the area and area ratio of aneurysm walls exposed to high/low WSS metrics determined by varying high threshold (HT) or low threshold (LT): (A) TAWSS LT , (B) TAWSS HT , (C) OSI HT , (D) WSSD LT , (E) WSSG HT , (F) GON HT .

TABLE 1
Comparison of WSS metrics computed by the global and local models.
*represents 0.01 < p < 0.05 and **represents p < 0.01.Frontiers in Bioengineering and Biotechnology frontiersin.orgstatisticalsignificance.The Spearman's rank test was conducted to evaluate the correlation between any two parameters.Multivariate logistic regression analyses were employed to derive models for discriminating the rupture status of aneurysms, with the performances of the models being evaluated by receiver operating characteristics (ROC) analysis.All the statistical analyses were implemented in SPSS 26 (IBM, Armonk, NY).

TABLE 3
Results of ROC analyses for individual morphological and hemodynamic parameters regarding their performances in discriminating the rupture status of aneurysms.The data were presented in form of AUC and sensitivity/specificity at an optimal threshold of each parameter.
Spearman's rank coefficients of the correlations between morphological/hemodynamic parameters of aneurysms.*indicates a p value < 0.01, and *indicates a p value in the range of 0.01 ≤ p < 0.05. *